Identification of vitamin D-related signature for predicting the clinical outcome and immunotherapy response in hepatocellular carcinoma

Hepatocellular carcinoma (HCC) is one of the most common cancers globally, seriously endangering people health. Vitamin D was significantly associated with tumor progression and patients’ prognosis. Integrative 10 machine learning algorithms were used to develop a Vitamin D-related signature (VRS) with one training cohort and 3 testing cohorts. The performance of VRS in predicting the immunology response was verified using several predicting approaches. The optimal VRS was constructed by stepCox + superPC algorithm. VRS acted as a risk factor for HCC patients. HCC patients with high-risk score had a poor clinical outcome and the AUCs of 1-, 3-, and 5-year ROC were 0.786, 0.755, and 0.786, respectively. A higher level of CD8 + cytotoxic T cells and B cells was obtained in HCC patients with low-risk score. There is higher PD1&CTLA4 immunophenoscore and TMB score in low-risk score in HCC patients. Lower TIDE score and tumor escape score was found in HCC cases with low-risk score. The IC50 value of camptothecin, docetaxel, crizotinib, dasatinib, and erlotinib was lower in HCC cases with high-risk score. HCC patients with high-risk score had a higher score of cancer-related hallmarks, including angiogenesis, glycolysis, and NOTCH signaling. Our study proposed a novel VRS for HCC, which served as an indicator for predicting clinical outcome and immunotherapy responses in HCC.


Introduction
Hepatocellular carcinoma (HCC) is one of the most common cancers globally, seriously endangering people health. [1]There are over 900,000 new diagnosed HCC cases and over 830,000 HCC-related deaths around the world in 2022. [2]Many HCC patients lose the opportunity for surgery, because they are usually diagnosed at advanced stages. [3]Less than 1/5 of HCC patients survive more than 5 years. [3]High heterogeneity, invasiveness and metastasis are the main reasons leading to the poor clinical outcome of HCC patients. [4]Limited biomarkers could be used for predicting drug sensitivity and clinical outcome clinically.
Vitamin D is a fat-soluble vitamin, [5] which could be obtained from diet or from sunlight-dependent endogenous synthesis of the lower epidermis. [6]Study showed that vitamin D could help control tumors and prolong the lives of cancer cases. [7]oreover, lack of vitamin D is correlated with higher cancer incidence rate. [8]Vitamin D suppressed tumor growth and involved in the development of cancer. [9]Leyssens et al found that vitamin D was involved in the proliferation and differentiation in CRC. [10]Vitamin D metabolism-related gene CYP24A1 was overexpressed in some types of cancer and correlated with tumor growth. [11,12]Considering the vital role of vitamin D in cancer, it is necessary to explore the role of vitamin D-related genes in HCC.
Based on bulk RNA-seq data, our study developed a vitamin D-related signature (VRS) for HCC cases.The correlation between VRS and drug sensitivity of HCC patients was evaluated by several predicting scores.

Data acquisition
Bulk RNA-seq datasets of HCC were obtained from the Cancer Genome Atlas (TCGA) (n = 329), ICGC (n = 228), GSE14520 (n = 218) and GSE72094 (n = 90) datasets.Metastatic HCC should be excluded from our study.IMvigor210 (n = 298) and GSE91061 (n = 89) dataset were utilized to explore the linked between the VRS and immunotherapy response.A total of 194 vitamin D-related gene lists were obtained from 2 literature pieces. [5,13]Among them, 173 vitamin D-related genes were identified from TCGA.

Integrative machine learning algorithms constructed an optimal VRS
DEGs were identified by "limma" package with the cutoff of |LogFC| ≥ 1.5.Risk factors were identified by univariate Cox analysis.An integrative analysis procedure with 10 machine learning algorithms was used to develop a predicting VRS for HCC.The candidate genes in VRS and corresponding coefficient were firstly obtained in TCGA cohort.The VRS score and C-index of HCC cases were then calculated in all cohorts.The prognostic signature with the highest average C-index was regarded as the optimal VRS.

Assessment of VRS
Based on the best cutoff determining by the survminer R package, HCC cases were separated into groups.We then drew ROC curve and C-index using the "survival" R package.The potential prognostic biomarkers were obtained after performing Cox analyses (Univariate and multivariate).The construction of predictive nomogram was determined by "nomogramEx" package using VRS score and other clinical parameters.To evaluate the possibility of VRS for clinical application, we then draw a decision curve analysis (DCA) curve.

Immune infiltration analysis
ESTIMATE analysis was performed to calculate ESTIMATEScore of HCC patients. [14]Immunedeconv (integrating 6 state-of-theart algorithms) was utilized to explore the correlation between VRS score and immune cell. [15]The level of immune cells, immune-related activities or functions scores were analyzed with ssGSEA method.The GSVA package was used to calculate the score of the gene set of "h.all.v7.4.symbols.gmt."

Drug sensitivity analysis
Several predicting scores (TMB score, immunophenoscore, tumor escape score and TIDE score) were utilized to explore the value of VRS in immunotherapy benefit predicting of HCC cases.There is less likelihood of immune escape and better immunotherapy benefit in patients with lower TIDE score and higher TMB score.Bu using the "oncoPredict" R package, we then obtained The IC50 of drugs in each HCC case

Evaluation of the performance of VRS
Univariate and multivariate Cox regression analysis suggested VRS-based risk score as an independent risk factor for the clinical outcome of HCC patients in TCGA, ICGC, GSE14520 and GSE76427 cohort (Fig. 3A-B).The C-index of VRS-based risk score was higher than that of age, gender and clinical stage (Fig. 3C-F).We also randomly collected 31 prognostic signatures that have constructed for HCC (Supplementary Table 1, http:// links.lww.com/MD/M326) and calculated their C-index.As a result, the C-index of VRS was higher than that most of these 31 prognostic signatures, excepted ChenY signature (Fig. 3G).To further evaluate the mortality risk at 1, 3 and 5 years of HCC patients, we then constructed a predictive nomogram (Fig. 3H).The calibration curves demonstrated the relative well conformity between speculated outcomes and observed outcomes (Fig. 3I).

VRS-based distinct ecosystem in HCC
Significant correlation was obtained between risk score and immune cells (Fig. 4A), including B cells, CD8 + T cells and dendritic cells (Fig. 4B-D).Moreover, a higher abundance of immunoactivated cells was obtained in HCC patients with low-risk score, including B cells, CD8 + T cells and NK cells (Fig. 4E).HCC patients with low-risk score had a higher stromal score, immune score and ESTIMAE score (Fig. 4F, all P < .05).We also found a higher cytolytic score, inflammation-promoting score and T cell co-stimulation score in HCC patients with low-risk score (Fig. 4G).Further analysis revealed that the level of most of the HLA-related genes was higher in HCC patients with low-risk group (Fig. 4H).

VRS acted as an indicator for the drug sensitivity in HCC
Several approaches were used to evaluate the predictive value of VRS-based risk score in immunotherapy response.The result  dysfunction score, T cell exclusion score (Fig. 5C, all P < .05).A lower immune escape score was obtained in HCC patients with low-risk (Fig. 5D, P = .045).Immune checkpoints played a vital role in inhibiting immune responses and promoted selftolerance in cancer.In our study, we found that the level of most of immune checkpoints was lower in HCC patients with low-risk score (Fig. 5E, all P < .05).This evidence showed that HCC with low-risk score may have a better sensitivity to immunotherapy.We also verified our results using 2 immunotherapy cohort.In GSE91061 cohort, the risk score in non-responders was higher and patients with high-risk score had a poor OS rate (Fig. 5F, P < .05).Moreover, the response rate in high-risk score group was significant lower (Fig. 5F).Similar results were obtained in IMvigor210 cohort (Fig. 5G).We also explored the IC50 value of common drugs for chemotherapy and targeted therapy in HCC.The result suggested a lower IC50 value of chemotherapy drug (5-fluorouracil, camptothecin, docetaxel, gemcitabine, Oxaliplatin, and paclitaxel) and targeted therapy drug (afatinib, crizotinib, dasatinib, erlotinib, and lapatinib) in HCC patients with high-risk score (Fig. 6A-B).Thus, HCC with high-risk score may be better sensitivity to chemotherapy and targeted therapy.

The distinct difference in cancer-related hallmarks in different VRS-based risk score group
We then performed gene set enrichment analysis to explore the biological process in HCC patients with different VRS score.As shown in Figure 7, HCC patients with high-risk score had a higher score of angiogenesis, coagulation, DNA repair, glycolysis, hypoxia, IL2_STAT5 signaling, MTORC1 signaling, and P53 pathway, indicating that the activation of these biological processes may play the vital role in HCC tumorigenesis and progression.

Discussion
In our study, a total of 10 integrative machine learning algorithms was used to construct a VRS, which was verified by using 3 validation cohort.The optimal VRS was constructed by stepCox + superPC algorithm.VRS acted as a risk factor for HCC patients.HCC patients with high-risk score had a poor clinical outcome and the AUCs of 1-, 3-, and 5-year ROC were 0.786, 0.755, and 0.786, respectively.Moreover, we also found that VRS acted as an indicator for predicting drug sensitivity in HCC.
In our study, we constructed a signature consisting of 12 VRGs which was prognostic biomarker of HCC, including CYP24A1, NCOA7, TGFB1, GSR, IGFBP3, IGFBP2, VGF, AGAP2, DENND6B, LRRC8A, BCL6, FCER2.TGFB1 was involved in the aggravation of HCC malignant behaviors. [16]GFBP3 inhibit cell growth and insulin-like growth factor signaling in HCC. [17]LRRC8A served as a central mediator and accelerated colon cancer metastasis by regulating PIP5K1B/PIP2 pathway. [18]NCOA7 could regulates growth and metastasis of clear cell renal cell carcinoma via MAPK/ ERK Signaling Pathway. [19]IGFBP3 promotes resistance to Olaparib via modulating EGFR signaling in advanced prostate cancer. [20]ne of the best therapy options for the patients that could not receive operation is immunotherapy. [21]However, there are limited sensitive markers for the prediction of immunotherapy benefits though many immune checkpoint-based drugs have been used clinically.In our study, the value of VRS in immunotherapy benefit prediction of patients.There is a higher immune escape likelihood and less immunotherapy benefits in patients with higher TIDE score. [22]Higher TMB score suggested betterimmunotherapy benefit. [23]We found a higher TMB score and lower TIDE score and tumor escape score in HCC patients with low-risk score.Immunophenotype score, primarily developed from TCGA RNA-seq data, was designed to predict patient responses to immune checkpoint inhibitor treatments. [24]It is confirmed that low-risk populations respond more effectively to immunotherapy, as evidenced by the higher immunophenotype score, which was consistent with our findings.Thus, low VRS-based risk score may indicate a batter response to immunotherapy.
We then performed gene set enrichment analysis to explore the biological process in HCC patients with different VRS score.The result suggested HCC patients with high-risk score had a higher score of angiogenesis, coagulation, glycolysis, hypoxia, MTORC1 signaling, and P53 pathway, indicating that the activation of these biological processes may play the vital role in HCC tumorigenesis and progression.Angiogenesis played a vital role in the developed HCC. [25]Hypoxia was also involved in innate immunity and tumor progression in HCC. [26]lycolysis shows significant correlation with the prognosis of HCC. [27]

Conclusion
Our study proposed a novel VRS for HCC, which served as an indicator for predicting clinical outcome and immunotherapy response in HCC.

Figure 1 .
Figure 1.Identification of potential prognostic biomarkers of vitamin D-related genes in HCC.(A) Heatmap showing the differentially expressed genes in HCC.(B) The overlap of vitamin D-related genes and differentially expressed genes in HCC.(C) The potential biomarkers of vitamin D-related genes in HCC.HCC = hepatocellular carcinoma.

Figure 2 .
Figure 2. Identification and evaluation of a Vitamin D-related prognostic signature (VRS).(A) The C-index of each prognostic model constructed by 10 machine learning algorithms in training and testing cohort.The survival curve and corresponding ROC curve of HCC with high and low-risk score in TCGA (B), ICGC (C), GSE14520 (D) and GSE76427 (E) cohort.HCC = hepatocellular carcinoma, TCGA = the Cancer Genome Atlas.

Figure 3 .
Figure 3. Evaluation of the performance of Vitamin D-related signature (VRS).Univariate (A) and multivariate (B) Cox regression analysis identified risk factor for HCC.(C-F) C-index curve evaluated the discrimination of VRS in predicting the overall survival rate of HCC patients in training and testing cohort.(G) The C-index of risk models had been developed for HCC.(H) A predictive nomogram constructing with risk score, age, gender and clinical stage.(I) Calibration plots demonstrated that the actual 1-yr, 3-yr, and 5-yr survival times were highly consistent with the predicted survival times.HCC = hepatocellular carcinoma.

Figure 4 .
Figure 4. Immune microenvironment landscape in HCC patients with different risk score.(A) The correlation between risk score and immune cell infiltration based on several state-of-the-art algorithms.Risk score showed negative correlation with the abundance of B cells (B), CD8 + T cells (C) and dendritic cells (D).(E) The level of immune cells in HCC patients with different risk score based on ssGSEA analysis.(F) Low-risk score indicated a higher stromal score, immune score, and ESTIMAE score.(G-H) The score of immune-related functions and HLA-related genes in HCC patients with different risk scores.*P < .05,**P < .01,***P < .001.HCC = hepatocellular carcinoma.

Figure 5 .
Figure 5. Vitamin D-related signature as an indicator for immunotherapy response in HCC.(A-B) The immunophenoscore and TMB score in HCC patients with high-and low-risk score.(C-D) The TIDE, T cell dysfunction score, T cell exclusion score and immune escape score in HCC patients with high-and lowrisk score.(E) The level of common immune checkpoints in HCC patients with high-and low-risk score.The overall rate and immunotherapy response rate in patients with high-and low-risk score in GSE91061 (F) and IMvigor210 (G) cohort.*P < .05,**P < .01,***P < .001.HCC = hepatocellular carcinoma.

Figure 6 .
Figure 6.The IC50 value of common drugs for chemotherapy and targeted therapy.High-risk score indicated a lower IC50 value of common drugs for chemotherapy (A) and target therapy (B) in HCC patients with high-risk score.HCC = hepatocellular carcinoma.